Energetics of metal slabs and clusters: the rectangle-box model 
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An expansion of energy characteristics of wide thin slab of thickness L in power of 1/L is con- 
structed using the free-electron approximation and the model of a potential well of finite depth. 
Accuracy of results in each order of the expansion is analyzed. Size dependences of the work func- 
tion and electronic elastic force for Au and Na slabs are calculated. It is concluded that the work 
function of low- dimensional metal structure is always smaller that of semi-infinite metal sample. 

A mechanism for the Coulomb instability of charged metal clusters, different from Rayleigh's 
one, is discussed. The two-component model of a metallic cluster yields the different critical sizes 
depending on a kind of charging particles (electrons or ions). For the cuboid clusters, the electronic 
spectrum quantization is taken into account. The calculated critical sizes of Ag^~ and Au^~ clusters 
are in a good agreement with experimental data. A qualitative explanation is suggested for the 
Coulomb explosion of positively charged NaJ^^ clusters at 3< n <5. 

PACS numbers: 73.21. La, 73.30.+y, 79.60.Dp, 68.65.La, 36.40.Qv 
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I. INTRODUCTION 

Clusters constitute a bridge between atomic, molecular 
and surface physics. Numerous investigations of physical 
properties of the low-dimensional systems are stipulated 
by their promising application in the nanotechnology. 

Recently 0, the work function of atomically uni- 
form Ag films grown on Fe(lOO) was measured as a 
function of film thickness. The maxima of the work 
function magnitude correspond to the "magic" thick- 
nesses. Scanning-tunneling microscopy observation of Pb 
nanocrystals grown on Cu(lll) indicates that in the equi- 
librium distribution of the island heights, some heights 
appear much more frequently than other ones ^ . An ap- 
pearance of these 'magic" island heights on the Cu fiat 
surface was studied by self-consistent electronic structure 
calculations 

Point contacts of gold bodies are investigated experi- 
mentally in Refs. 0,13] during elongation of the contacts 
to the rupture. It is shown that oscillations of elastic 
constants appear simultaneously with an abrupt change 
in conductance. A dimensionality of a contact varies dur- 
ing the process of stretching. Indeed, at the moment of 
formation of a contact, the contact region can be repre- 
sented as a slab inserted between electrodes, whereas at 
the moment of its rupture, the contact region becomes 
a wire. Thus, in the experiment we observe a transition 
from 2D (or OD) to ID open electron system. 

Analytical approaches to the determination of the den- 
sity of states and the Fermi energy for metal slabs are 
proposed in Refs. 13,13] on the basis of the free electron 
model. Quantization of a contact potential difference for 



a slab was described in Ref. (8| within the framework 
of the model of hard walls. Jumps of the tensile force 
in a point contact were explained in Refs. 0, 0, ITU . 
However, the work function cannot be determined in this 
simplest model (see review [T^). 

Probably, the electron work function of a slab was cal- 
culated for the first time by Schulte ^3 ■ The work func- 
tion ma gnit ude showed oscillations near its average value. 
In Ref. (Tj], the above approaches have been criticised. 
Detailed computations J,5, 16, 17, 18, 19, 20, 21] (includ- 
ing ab initio calculations) performed to date do not yield 
an unequivocal conclusion about size dependence of the 
work function of isolated slabs and wires. In addition, 
amplitudes of the work function oscillations are larger 
than in experiment. 

Since the work of Sattler et al. [22], mass- 
spectrometric investigations of the charging effects in 
cluster beams have clearly demonstrated the size- 
dependent Coulomb instability of clusters composed of 
a countable number of atoms [23, 123 . |25| . 

Rayleigh's theory predicts an instability of a charged 
liquid sphere of R radius, when the Hartree energy ex- 
ceeds twice the surface energy. The critical charge is 
defined by the expression 



= ±\/l67ri?3r. 



(1) 
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where r is the surface tension (or stress). Recently, this 
criterium of stability has been confirmed in the experi- 
ment 26] for microdroplets of ethylene glycol. 

Eq. 0J (i.e. Rayleigh's criterium) does not deter- 
mine, which type of the particles charge the cluster. A 
metallic droplet can contain either an excess number of 
electrons — \Qb\/& or ions AiV^ — jQaj/Ze, where 

Z is the valence and e is the elementary positive charge. 
Therefore, such a problem should be considered using 
the two-component model of a cluster in which electrons 
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and ions are interpreted on equal footing [23, ■ A sign 
of the excess charge resuhs in different size dependence, 
ATVe.i cx i? or ATVr cx R^'"^. 

In the present work, we developed an analytical the- 
ory of size-dependent energy and force characteristics for 
metal slabs using an elementary one-particle approach. 
This simple model makes it possible to calculate the os- 
cillatory size dependence of the work function and elastic 
force. Thermal effects are not considered. An assump- 
tion of the ideal plastic strain (when volume of the slab 
remains constant during stretching) allows for the com- 
parison of theoretical results and experimental ones 

In the framework of the two-component model, the 
size-dependent Coulomb instability of charged metallic 
clusters was described. A spectrum quantization is taken 
into account for a cluster having the shape of paral- 
lelepiped. The model makes it possible to study a physi- 
cal origin of the instability. Theoretical results on critical 
sizes of different shaped charged Au, Ag and Na clus- 
ters are in agreement with experimental ones (see Refs. 
[23LI25II and referencee therein). 



II. SLABS: QUASICLASSICAL 
APPROXIMATION 



we use the well-known expression 

where j, s, i are the natural numbers. 

It is convenient to use dimensionless variables by choos- 
ing Uo as a unit of energy and as a unit of length. 
We introduce the following notation: 

— kf)L^/2TT, ly = /co^y/STT, I ~ k^Lz/iT. 

Note that energy can be interpreted as the square of the 
state vector in the ^— space, = + + with 
< 1 • Eq- 101 gets a form 

2 

^^i = arcsin^i -I- i. (3) 

TT 

Not only solutions of Eq. lO but also their number are 
fully determined by the thickness Z, namely, ?f = W + 1; 
where [a] denotes the integer part of a. 

B. Density of states 



A. Formulation of problem 

We consider a thin slab with the thickness, — Apo 
(Apo is the Fermi wavelength of electrons in a semi- 
infinite metal) . The thickness is much smaller than other 
dimensions, ^ L^,Ly. Thus, the discreteness of the 
electron momentum components Px and py can be ig- 
nored. For typical values of the electron concentration in 
metals, we have Apo — 0.5 nm. 

As a first approximation, the profile of the one-electron 
effective potential in the slab can be represented as a 
rectangular potential well of constant depth Uq <Q with 
dimensions Lx , Ly, and L^. A solution of the three- 
dimensional Schrodinger equation for a quantum box is 
simple. An equation can be decoupled into three one- 
dimensional equations. Therefore it is characterized by 
a set of electron wave numbers /cxj — Sttj/Lx, fcys — 
27rs/Ly, and fczi, which are roots of the equation 

^zi-^z = — 2 arcsin(fczi/fco) + ttz, (2) 



where feg = •\/2m|C/o| and m is the electron mass. The 
numbers j, s — 0, ±1, ±2, ±3, ... and i = 1, 2, 3, ... 

The set of wave numbers determines the electron en- 
ergy: 

where p is the number of the electron state (the states are 
numbered in order of increasing of the electron energy) . 
For a cuboid with hard walls and dimensions Lx, Ly, L^, 



Let us estimate the interval A^ between neighboring 
values of ^z- It follows from Eq. that, for sufficiently 
large values of /, 

Ae « l/l. (4) 

Distances between two consecutive values of ^x and are 
small, namely, A^x = Cxj+i — ■Cxj = l/^x and A^y = 1/ly. 
For the relative values of the slab dimensions assumed 
by us, it can be easily found that A<^ ^ A^x, Ai^y. We 
see that the electron states {^xj, Cys, £,[} form a system of 
parallel planes ^z = 6 in ^— space and that the density 
of states on all these planes is the same and equal to 

a = 2/(AexACy) = 21 Jy. (5) 

The factor 2 takes into account two possible values of the 
electron spin polarization. 

Electrons occupy the states, beginning from the point 
{0,0,^1}, in ascending order of ^p, i.e., of the state en- 
ergy. Therefore, it appears that all occupied states lie in 
the ^-space domain bounded by the plane ^z = Ci and 
the hemisphere of radius ^p = ^/EF/\Uo\, where E^p > 
is the Fermi energy equal to the maximum energy of oc- 
cupied states. 

The occupied states are distributed with density a over 
disks formed by intersections of the Fermi hemisphere 
with the planes ^z — i = 1,2, ...,zp. The area of the 
disk is Si = 7r(i^p — The number of occupied states 
coincides with the number of free electrons in the slab, 

7Ve = af]5i = 2Uyf]7r(^|-C?), (6) 

i=l i=l 
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where ip is the number of roots of Eq. 

The number of occupied states per unit volume is 



3l 




(7) 



where we used the notation ep = £'p/|J7o| — £,f- 

By definition, the density of states p{E) is the num- 
ber of states per unit energy interval near the energy E 
and per unit volume of the metal. In order to find this 
quantity, we write Eq. in the form 



i=l / 



(8) 



One can interpret v as a number of states (per unit vol- 
ume) whose energies do not exceed e. In Eq. ((SJ is the 
index of the greatest of the roots of Eq. Q satisfying 
the condition ^? < e. 

We find from Eq. © that 



z = /^i H — arcsin^j. 

TT 

Substituting here by ^, we let ^ to take any value in 
the limits from to 1 and form an integer- valued in- 
creasing function i(^) such that at the points £,—£,[ the 
value of the function is increased by one and in the inter- 
vals between these points the function docs not change. 
Substituting ^ = we obtain 



ly/^ H — arcsin\/e 

TT 



(9) 



Square brackets indicate an integer part. 

By differentiating l.h.s. of Eq. © with respect to 
under the condition = const, we find p{e). Using Eq. 

and working backward through the normalizations, 
we obtain 



PiE) 



1 dNe 
V~dE 



where V — L-xL^Ly 




TT71 > (10) 



C. Size dependence of the work function 
From Eq. JT)), we obtain 



Using Eq. 0, we find that 



i=l 



IF 



l^/eF H — arcsin ^Jep 

TT 



(11) 



(12) 



In what follows, we assume that the electron density 
in the slab does not depend on its size and is 



V 



1,3 

47r3 



(13) 



If the depth of the well is fixed, then Eq. (|13|) implies 
that V = const. The thickness dependence of the Fermi 
energy ep(Z) can be found by solving the set of equations 
(lllll and p2(l under the additional condition v = const. 



Substituting the expression 



= fcpo/BTT^ into Eq. 



()lr{|l (fcpo is the Fermi wave number for the semi-infinite 
metal), we find 



3 

FO' 



(14) 



where ^po = /cpo/fcoi i-e. v is equal to doubled volume 
of the Fermi hemisphere in ^— space in the limiting case 
I — > oo. 

Setting ip — const, from Eq. Hll|l we obtain 



(15) 



In order to see how the roots of Eq. © change with 
varying I, we differentiate both parts of this equation 
and find that 




— £ = 

dn 1 + 



2^f 



< 0. 



(16) 



oo. I.e., as 



Here, the equality valid only in the limit I 
£i ^ for all i. 

The disks £z — £i are lowered with increasing I. The 
lowering rate decreases gradually, so that the lower disks 
move more slowly that the higher ones. Accordingly, the 
distance between the disks decreases and their number 
ip grows. It is seen from Eq. H12|l that this number in- 
creases by 1 each time the equality ep = Cip+i satisfied. 
The process of disk lowering is accompanied by the "pul- 
sation" of the Fermi hemisphere. Its radius £f = ^/£f 
alternately increases (as dep/dl > 0), and decreases ( as 
dep /dl < 0), having the average tendency to decrease. At 
S.F = the derivative dep/dl is discontinuous. The 

value of the jump decreases with increasing I. 

The minimum value of corresponds to the thickness 
equal to the atom diameter. Let us estimate the mini- 
mum value of I. We use I — L2y^2ni\Uo\/{TTh), where 
— 0.5 nm (it means that we have a single layer of atoms), 
and 



|C/o| =£;po + Wco, Ei 



FO 



2m 



(17) 



The work function for a semi-indefinite metal equals 
to 2.25 and 4.25 eV for Cs and Al, respectively |29l|. As a 
result, we have 1.6 < Zniin < 3.5. We assume henceforth 
the value l/l to be small and apply an expansion in terms 
of l/l for calculating of the Fermi energy. 
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Now we introduce the notation a = l/l. The S,i{a) 
dependence is imphcitly determined by Eq. which 
can be written as 



6 2 

— = arcsm^i + i. 

a TT 



(18) 



We look for the roots in the form of the expansion 

6 = 6 l„o + C( L.o « + L.o o? + i^r l„o o? + ... 

(19) 

Keeping the terms to the order of , we obtain the Fermi 
energy up to the order . It is seen from Eq. H18() that 



(20) 



Differentiating both sides of Eq. (|18|l with respect to a 
and multiplying the result by a, we obtain 



1 



a 



Setting a — Q, we find 

In a similar way, we obtain 



(21) 



(22) 



4? 

TT 



24i 



Substituting the obtained expressions into Eq. (|19|l . we 

find 



6 = ~ — + + ©(a"*). 

TT TT^ 



(23) 



Now we evaluate the Fermi energy to the first order 
in a. For this purpose, it suffices to keep the first two 
terms in Eq. (|23|l when substituting it into Eq. (|ll(l . 
Indeed, the order of magnitude of the error 5S^i does not 
exceed ipa^ in this case. The error of is 2^i(5^i, and 
its order of magnitude also does not exceed iyo? , since 

< 1. The error of the sum X^lLi is smaller than 
ipO;^, and the order of magnitude of the error of the whole 
expression (|11|) does not exceed i^ct' ■ Since zp — ^Je^-I 
this estimation of the Fermi energy is correct to first order 
in a. 

After the above-mentioned substitution into Eq. (|ll|l . 
we have 



£F 



V 1 

27r iTto. 



F „3 



a'' + 0(a'). (24) 



We divide the range of variation of a into intervals 
(cki+i,Q;i), i = 2,3..., so that = l/Zj we have ip — i 
inside these intervals. The values a > 0.3, which cor- 
respond to I < Imin, are nonphysical. Let us find the 
boundaries of the intervals. 

From Eqs. H23|l . (|24|l and condition ep — we can 
obtain an equation for Oi+i: 



8?^ 
3^ 



*i+i 



"3~ 



3f 



(25) 



According to Descartes' rule of signs, Eq. H25() has two 
real positive roots. In zeroth approximation, one of them 
is ^ ^/spo/i and the other is of a higher order of small- 
ness. We are interested in the first root of Eq. H25|l . be- 
cause boundaries of the interval (ai+i,ai), with constant 
ip inside, are determined by this root. The boundaries 
are 



(26) 



The minus sign corresponds to p = i + 1, and the plus 
sign to p — i. Thus, the width of the interval (ai+i, a\) 
decreases as l/i^ with increasing i (or I), since cti — ai+i = 

Unexpectedly, the examination of function £p{a) in- 
dicates the insufficiency of its approximation by the ex- 
pression (|24|l . This expression is stair- like function, i.e., 
it is constant within intervals (ai+i, ai), where ip = i. It 
can be explained by the fact that ai+i and ai differ in 
the second order of smallness, and after substituting into 
they give the same results (of course, terms of the 
order more than first must be omitted in the resultant 
expression) . A numerical calculation of ep by the formula 
(|24|l leads to the error that is associated with an incorrect 
consideration for terms of the order more than first and 
gives an appearance of the oscillatory dependence eF(^). 



In order to take the terms 



on the right-hand side 



of Eq. 1)24(1 into account, we must include the following 
terms, which were earlier neglected: 



2i 



-a 



4,- 2 



(27) 



In this case, the £f(0 dependence is represented by a 
concave curve in each interval (/i,Zi+i), i = 2,3, ... (/i = 

At the points I = h, the derivative dep/dl is discontin- 
uous and its jump is — 2eFo'^^^/i^. To the left from this 
point, the function grows, and to the right from it, the 
function decreases. The jump in the derivative results in 
the appearance of cusps in the plot. The sharpness of 
the cusps decreases with increasing i. 

For large values of I, Eq. (|12|l can be written as 



ip 



0{a^). 



(28) 



Substitute Eq. (^5)) into Eq. and using conven- 

tional units, we obtain 



2 V 2m \^ 37r\/ \Uo\ ) ^ ' 

The expression in the brackets is positive, i.e., asymptot- 
ically, we always have Ep > EpQ. 

In this model, the work function is defined as: 



= -Un - Ep 



(30) 



and it is easy to see that Wo < W^o- A role of the size 
dependence of the bottom of the potential well U{L.^) will 
be discussed later. 
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D. Deformation force 

In order to calculate force characteristics, we must find 
the size dependent electron kinetic energy. We denote the 
total kinetic energy of the electrons by e = K/\Uo\. 

As noted above, the one-electron kinetic energy e-p is 
numerically equal to the square of the radius-vector of 
the point in ^— space. The contribution from the cor- 
responding element dS of the disk to the total kinetic 
energy is de = EpCrdS, where the density of states a is 
defined by Eq. ISJ- Next, we must integrate over the 
disk area and sum the contributions from all disks. 

We introduce the notation p = (^^j -I- ^y^)^/^. The 
maximum value of p in the i-th disk is equal to the disk 
radius pi = (ep — We have 



(31) 



Performing the summation, we find 

e = TTlJy —T- — IT + —i-pa - — ipCk 
^ lYO^ Stt 45 

^ ■ .2 2 1-4 4 32 .,- ,-\ , 

— zpa-^^iptt - -i-po. + T^ipa ) ]+U[a). 



27r 



6 



457r 



(32) 



The asymptotic form of this expression in conventional 
units is 



K = -NrE- 



FO 




1 - 




Let us discuss the origin of the different terms in Eq. 
(ESI. 

The first term in Eq. I|33|l is the kinetic energy in the 
case, when the slab thickness is comparable to other 
dimensions. The distance between the disks in ^— space 
in this case is so small that the summation in Eq. H31|) 
can be replaced by integration. 

The second term in the brackets in Eq. (|33|l appears for 
a finite well depth. This correction is rather important 
and makes a contribution of about 50%. In contrast to 
the case of an infinite well, electron localization in a well 
of finite depth is not strict, therefore, the kinetic energy 
in the latter case is smaller. 

In order to compare our results with the experimental 
data (4,, , we find the oscillating electron contribution 
to the elastic force under the conditions of ideal plastic 
strains, i.e., in the case where total volume of the slab is 
conserved: 



/z 



-{ds/dh)y. 



This part of the force has no relation to the phases of 
stretching that are accompanied by a change in volume. 
It rather determines the variation in slab elastic proper- 
ties as the slab thickness is varied. This force depends 
on the number of particles in the slab, therefore, it is 
convenient to consider the force normalized by iVo, 



N, 



2™ 



3 LI 
57r3 
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(34) 



III. CHARGED CLUSTERS 

A. Two-component model: asymptotic expressions 

Let us consider now a neutral cluster, which contains 
Nc/Z = Ni = N atoms. Total energy of a cluster, 
charged by AiVc ^ electrons, can be written as [30|| : 



Ne+AN„ 



2C 



(35) 



where /Ze is the electron chemical potential. A cluster 
can retain the excess electrons ANc only when its energy 
in this state is lower than the energy in the state with 
Nc + ANe — 1 electrons. By our definition, the number of 
electrons in a cluster is critical, if the number of electrons 



AiV*, for which the reaction 



(an:)- 



(an;;-!)- 



is reversible, and the ionization potential of the cluster, 
IP* > 0, tends to zero 



A^(AiV;) = ^N,+ANJ 



N„-|-AN* 



IP* 



0. (36) 



We note that the addition of one more surplus electrons 
to AN* is possible only in a metastable state, because 
the affinity of this electron is 



EA* = El 



E^ 



hANj-l-l 
„2 



— (2AAr,-t-I)<0. (37) 



In any case, the relation 



IP* - EA* = — 

is valid. When AA^o > AA"*, the cluster is overcharged. 
The surplus electrons are separated from the free states 
by a potential barrier and they can be in the metastable 
state for some time. The lifetime of each electron is deter- 
mined by specific conditions in the nonequilibrium sys- 
tem. 
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Using Eqs. (|36(l and (|35|) . we obtain for the critical 
excess electron charge 



an: 



WcqC - flcl 



(38) 



where Woo — — McO is the electron work function for a flat 
surface, /Ze = ^ieO + fJ'ci/R, and fi^i/R is the first curva- 
ture correction term in the electron chemical potential of 
a metal sphere of i? = iV^/'^ro radius with the notation 
To for the mean ion spacing. 

Let us consider now a positively charged cluster of 
metal, which contains — ZN{ electrons and A^i + AiVi 
ions. This situation is similar to that one when the 
droplet with A^i ions "contains" AiVo < lacking elec- 
trons. In this connection, AiVo has to be divisible by 
Z. 

The energy of charged cluster i^Ni+ANi is related to 
the energy of a neutral cluster in the following way 



Ni-t-ANi 



i+eZAN;) 
2C 



(39) 



As in Eq. H35|) the most essential size dependence is 
due to the self-repulsion of the surface charge -+eZlS.Ni. 
Actually, the ions are not mobile and the redistribution 
in the electronic subsystem mimics the distribution of an 
excess positive charge. 

The change in the total energy, associated with the 
detachment of AA^j— th ion, is 



AE{AN;) = i?N.+ANi-l 



-^Ni+ANi 



= -Mi 



2C 



(2A7Vi - 1). (40) 



A cluster with charge +eZAA^i can exist in equilibrium 
only if AE{ANi) > 0. So the number of ions AN* in a 
cluster is critical, if the reaction 



-"^Ni+AN* 



" ^"Ni+AN*-1 +^'-'1 



is reversible. In this case we have 



AiV* 



[Zef 



1 

2' 



(41) 



where Wio = — /iio is the ion work function for a plane 
surface. For a sphere of radius R = (iV; + AN^Y^^rQ, 
using sum rules [23, we can write 



27 

/^ii = — 

no. 



(42) 



where 7 is the surface energy per unit area. For materials 
under the study ^oi ^ 1-9 eVxoo [2II3H . 

If ANi > AN^ , the cluster emits a surplus ion, pass- 
ing into the state with lower energy. This approach cor- 
responds to the consideration of the droplet as a two- 
component electron-ion system. 



The work function of single-charged ion can be ex- 
pressed by means of the Born cycle js^l using the ion- 
ization potential of single atom IP{1), cohesive energy 
£cohOj and work function for a flat surface Wcq'- 



W^n - e. 



cohO 



/P(l) - 1^00- 



(43) 



For Pb Ecoho = 1-5 eV, Wco = 4.0 eV, IP{1) = 7.4 eV, and 
we get Wio = 4.9 eV. When R = 12ao, the critical charge 
is equal to +2.7e. This estimation is in a good agreement 
with experimental data and with the results of more 
complicated self-consistent calculations p^ . 

Our approximation assumes that the shape of a clus- 
ter is unvariable during its charging. Eq. H38|l for AA^* 
and Eq. (|41|l for AA^j* take into account electron and 
ion emission and distinguish between them. It is due to 
the necessity to expend the energy on the embedding of 
a particle of this kind into the cluster and on the redis- 
tribution of its charge over the surface. Such a mech- 
anism of the Coulomb instability can be introduced as 
an alternative to Rayleigh's one. Estimates show that 
AA^R > AA^j* > AA^*, i.e. during the charging, rather 
the single-electron/ion emission occurs than Rayleighs in- 
stability. For small clusters, the energy quantization is 
important. 



B. Quantum spectra 

The shape of real charged clusters seldom looks like a 
sphere (e.g. see |33,0l)i therefore it is convenient to de- 
termine the electron spectrum using the parallelepiped 
model. The allowable levels form a discrete spectrum. 
Wave vector components are evaluated by solutions of 
Eq. (PI for each directions. In order to separate real 
levels from virtual ones we introduce the following con- 
dition: 



kp/ko < 1. 



(44) 



Using the perturbation theory, Eq. ^ can be reduced 
to the infinite well problem |3^ . Three relations of iden- 
tical form determine wave numbers, for example. 



"■XJ — n-xj 



Afcxj, C=|A^xj/fc^|«l, (45) 



where k^ = nj/L^ is the solution at ko — > 00. Sub- 
stituting Eqs. H45|l into Eq. jSJ, we obtain in the first 
approximation for a cube C = — 2/fcoL and the energy 
spectrum 



2mL'' 



.(j2 + s2 + z2)(l + 2C + 0(C')). 



(46) 



One more alternative expression can be derived from Eq. 
(2) under the condition (44): 



E^ 



2m 



ko 



kf)L 



if 



(47) 
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2 6 10 14 18 

Figure 1: Size dependence of the work function for metallic 
slabs. 



On the one hand, in a neutral cube number of elec- 
trons is given, on the other hand it is determined by sum 
2^p(5(-E — Ep) over all occupied levels taking into ac- 
count twofold spin degeneracy. Filling levels by electrons 
we find a highest occupied state, < 0, counted off 

from the vacuum level. Following the Koopmans' theo- 
rem the ionization potential of a cubiform cluster can be 
determined as 

IP = -E^^ + ^. (48) 



IV. RESULTS AND DISCUSSION 

We performed calculation for slabs of trivalent Al, 
monovalent Au and Na with electron concentration fie — 
3/47rrg, where rg — 2.07, 3.01 and 3.99 oq, respectively. 
The work function values for semi- infinite metals, which 
we used, are W^o = 4.25, 5.15 (or 4.3) and 2.7 eV 36|. 

Figure |21 shows the results of calculation of the thick- 
ness dependence of the work function for extended iso- 
lated slabs. The inequality Wo < Woo is satisfied for the 
whole thickness range. The values of the largest oscilla- 
tions of the work function are about 0.1 - 0.2 eV. This 
size dependence is in a general agreement with the ex- 
perimental results and self-consistent calculations for 
extended thin Al slabs [IQl and cylindrical Al and Na 
wires [19II21I . However, there is a disagreement with the 
results of [H El HI Eg . 

Comparing We{L) for various metals, it is easily seen 
that all differences are determined by the values of . For 
aluminum (with the smallest rg), the amplitude of oscil- 
lations of the work function Wc/Woo is the largest, the 
period AL is the smallest, and position Li of the cusps are 
displaced to the left. These features are well described by 
the approximate relations 1 — Wc/Wco ~ l/vsL, AL ~ rg, 
and Li irg {i is the number of subband), which follow 
from Eqs. ^ and 



Table L The energy of the first occupied state, £"0 < 0, in 
spherical clusters Naw calculated by Ekardt |3^ . 



N 


18 20 34 40 58 68 90 


-Eo [eV] 


5.10 5.15 5.41 5.41 5.59 5.53 5.76 


N 


92 106 132 138 168 186 198 


-Eo [cV] 


5.63 5.64 5.86 5.79 5.91 5.92 5.81 



In order to clarify the role of the thickness dependence 
of the bottom of the well (see Eq. (|17|l ). we consider the 
data presented in Table |l| These data are extracted from 
the results of the self-consistent calculations performed 
in [3^. In that work, the electron energy spectrum was 
calculated for a self-consistent spherical potential whose 
form was far from being rectangular. In the rectangle-box 
model the position of the first occupied level Eq varies 
in accordance with the position of flat bottom of the po- 
tential well. As the well width is large (i.e., as ^ 00), 
this level is "lowered" to the bottom. Hence, we can ob- 
tain quite reliable information about the size dependence 
of the depth of the rectangular well by determining the 
confinement behavior of the lowest level in the potential 
profile corresponding to the spherical cluster. 

This dependence is almost monotonic and asymptot- 
ically weak. Moreover, it does not compete with the 
thickness dependence of the Fermi energy in Eq. H30(l 
and gives a minor contribution to Eq. H29|l . Taking 
into account the dependence U{L) in Eq. (|30|l leads to 
strengthening of the inequality Wc < Wcq. 

When the slab is inserted into contact with the elec- 
trodes, the electron chemical potentials are equalized and 
the electronic system should be considered as an open 
system with We{L) = Wbq. The electrical neutrality of 
the slab cluster is broken, and the part SN^ > of the 
electron liquid passes out to the bath. As a result, a 
contact potential difference Scj) appears. 

In order to determine the contact potential difference, 
we consider energy cycles, in which electrons are trans- 
ferred first to infinity and then to the electrodes. By 
analogy with Eq. H36(l . we express the ionization poten- 
tial of the slab, having charge +e6Ne, as 

IP — En^-SN^-A — ^^N„-5N„ 

= W,A + ^iiSN, + Af-6N^) (49) 

and write the electron affinity of the bath for the charge 
— eA as EA — WcqA. Equating these two quantities, we 
obtain 

W,o-W,- ^i2dN, + A) = 0, (50) 

where C is the sample capacitance. 

We note that A can be infinitesimal, since an elec- 
tron can pass through the contact only partially (i.e., it 
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Figure 2: Size dependence of the oscillating component of 
elastic force F^/Nc for metallic slabs. 



can be detected on both sides of the geometrical contact 
with a nonzero probabihty). JTVo can be considered as a 
smoothly varying quantity. This situation is typical for 
one-electron devices [s^ . 

We also assume that C corresponds to the total ca- 
pacitance of the both contacts. The validity of this as- 
sumption depends on the sample geometry and electro- 
magnetic environment. This is not true for a spherical 
cluster in contact with electrodes, but it is valid for a 
cubiform cluster or a slab |39l |. Neglecting environment 
and setting C = eSNc/Scj), dN^ < iVe, and A ^ 0, we 
obtain from Eq. l(5n|) 

5^ = (Wco - Wo)/e. (51) 

Now the energy spectrum of the Ni = — 6Ne elec- 
trons that remain in the slab must be found for the rect- 
angular potential well of changed depth Uq — eS(j). The 
total kinetic energy Ki of the remaining electrons can 
be determined in the same way as for an isolated slab 
but with changed energy spectrum and number of elec- 
trons. For the oscillating part of the elastic force, we 
have Fzi = — {dVt/ dL^)^, where Vt — Ki + WcqA^i is the 
thermodynamic potential. 

Using the data from Fig. 1, we can also determine 
the contact potential difference 5(t). This potential dif- 
ference produces a negative shift of the well depth (for 
the thinnest slab, it attains 0.5 - 1 eV). This leads to 
a displacement of density of states to values correspond- 
ing to greater thicknesses. The stretched sample acts as 
an "electron pump" with respect to the electrodes, alter- 
nately ejecting electronic liquid and drawing it back. 

In Fig. 121 we show the part of the elastic force that 
is due to the quantization. As can be seen from this 
figure, the amplitude of the force oscillations depends 
strongly on Vg. For sodium, the oscillation amplitude is 
8 times smaller than that for Aluminum. The character 
of this dependence can be determined from Eq. (|34|l : 
(Fz/iVe)i ~ ir^. Our estimations showed that, for a slab 
in a contact, force oscillations are analogous in shape and 
amplitude. 




0.1 0.2 0.3 0.4 0.5 



Figure 3: Size dependences of the first ionization potential 
14811 of AIn clusters. Solid, dashed, and dot-dashed lines indi- 
cate the dependences determined according to solutions Eqs. 
@) 14611 ■ and 1471 1. respectively. Numbers at the top corre- 
spond to the numbers of atoms in the cluster. 



The first maximum of the oscillating part of the force 
F^/Ne for Au (corresponding to a slab thickness of one 
monolayer) is 0.2 nN, i.e., it is much smaller than the 
experimental value for a wire, which is equal to 1.5 nN j^. 
This difference can be explained by both the difference 
in the dimensionality of the electron gas in a wire and in 
a slab and the effect of current flowing through a contact 
estimated in Ref. jl^l and discussed in Ref. ^l|- Note 
that there are no experimental results for slabs. 

For the case of aluminum cluster, let us apply analyti- 
cal approach of the previous section. Size dependences of 
the ionization potential l|48|l calculated for the cubiform 
clusters AIn are shown in Fig. O Here we use spectra, 
which are determined by Eqs. lO, and (|T7jl . For 

the A'' range (10, 3000), calculations designate an essen- 
tial role of spectrum quantization even for large clusters. 
Magic numbers obtained are close to those found exper- 
imentally They are different in the case of single 
valence sodium clusters. 

Beginning from one hundred atoms, calculations of the 
spectrum based on approximate expressions (|46|l and 147|) 
give quite reasonable results. However, their inaccuracy 
results in the level hierarchy differing from that deter- 
mined from Eq. There is a difference between spec- 
tra calculated by Eqs. H46I) and (|47|l for N values near 
A^ = 58 (see insert in Fig. 3) 43]. 

At the next step we investigate the ionization potential 
of a cluster as a function of its shape. We assume that 
shape of a cuboid cluster varies from strongly flattened to 
elongated. Thus, we have a monatomic slab of thickness 
L at the beginning and a monatomic chain of length L 
at the end. During this evolution, the sample volume is 
supposed to be constant and equal to 4 nm''. 

We divide the range of the size variation into 10'^ in- 
tervals and find the spectrum for each of them using Eq. 
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Figure 4: Evolution in the size dependence of the first ioniza- 
tion potential 1481 1 of NaN cubiform cluster with a change in 
the cluster shape from a slab to a wire. Dashed lines repre- 
sent the bottom position of the potential well, ?7o, and work 
function for the flat surface, VKeO- 



^ . We exchange a capacitance of the parallelepiped by 
capacitance of equivalent spheroid in the ionization po- 
tential (|48|1 . The size dependence of the capacitance has 
a minimum for a sphere. In the limiting cases of the slab 
and the wire of monoatomic thickness capacitances are 
approximately twice and sevenfold larger, respectively. 

Figure 0] displays the behavior of the electron work 
function and the ionization potential of the isolated 
sodium samples of varying shapes. The inequality 
_^HO ^ Wf^Q is observed to be obeyed over the whole 
range of the considered dimensionalities. The size depen- 
dence of energy of the first occupied state, Eo{L), has 
minimum at the point corresponding to a cubiform clus- 
ter. As it can be seen, there are ranges of lengths, where 
IP < Wco and IP > Weo. The inequality IP < W^o 
is rather surprising. It is known from experiments that 
the work function Weo of alkali metals is approximately 
equal to the one half of IP of the atom |4J|. Therefore 
one would expect that the value of IP of a small solid 
sample belongs to the interval Wcq < IP < /P(atom) 
independently of the shape of sample surface. However, 
the competition between the size correction Wc{L) and 
the e^/2C term in expression (|48|l for the ionization po- 
tential can lead to the opposite inequality. 

Using mass-spectrometer p5| . minimal number of 
atoms N have been recently determined for which an ex- 
istence of stable charged clusters Au^~ (TV >27), Au^~ 
{N >58), and Ag^~ {N >27) is possible under the condi- 
tion of a particle retaining 2 or 3 surplus electrons. This 
problem is inverse to the one considered above. Here 
AiV* is a parameter and N is unknown. 

We note that even for particles containing more than 
a thousand of ions the critical charge does not exceed 
a few units. This interesting feature is a result of the 
strong Coulomb repulsion of the surplus charge spread 



over the surface of the particle. The situation is different 
for atomic and molecular ions, where electrons are not 
collectivized. 

For calculations, we use the following experimental val- 
ues of surface tension: t = 1134, 780, 191 erg/cm^ for 
Au, Ag, and Na, respectively. Here, we assume that the 
surface energy equals the surface tension, while in reality 
these values can be considerably different from [isj . 

Rayleigh's expression gives numbers of atoms in the 
critical clusters, which are 4-5 times less. These num- 
bers are N ^ 9 and 6, for Au^~ and Ag^~, respectively. 
In our approach the problem is reduced to the solution 
of the equation 



IP*{AN*,N) 



{an:,n) 

(2AiV; - 1) = 0. (52) 



2Ceff(iV) 



Effective capacitance Ccff = i? -I- ^ is used in order to 
explain experimental results for charged clusters. The 
additional small quantity 6 is caused by an increase of 
radius of the charging electron "cloud". The 5 value 
was introduced for calculations of the polarization by 
Snider and Sorbello l46[ and the ionization potential 
of clusters by Perdew |47l |. The averaged dependence 
(5(ro) = 1.617 + 0.199(ro/Zi/3 - 2.07) (in bohrs) is ob- 
tained using coordinates of the image plane for various 
crystallographic faces l2a|. which were computed in the 
stabilized jellium model pSj . 

Note that introduction of 6 to Eqs. (|^ and lO 
is not rigorous. This procedure accounts only for the 
Hartree contribution, S/R^, into correction term ~ l/i?^ 
in the energy expansion. Nevertheless, solution of 

Eq. ifCT) is responsible to the S value With dis- 

cussed modification for IP, Eq. H36() is adequate for ex- 
plaining consecutive fotoionization acts for large clusters 
AIn in the wide N range (2000, 32000) The use 

of Ccff = R + S makes the size monotonic component of 
IP{N) in Fig. 3 some weaker. 

The size dependence IP*{AN*, N) calculated from re- 
lationships (|52|l and (|38|l is shown in Fig. |31 Points of 
the diagrams with IP* = indicate TV values for the 
critical clusters. It can be seen from Fig. 5 that the qua- 
siclassical dependence (|38|l and Eq. (|52|l with inclusion of 
the level quantization leads both to the better agreement 
with experimental data than Rayleigh's formula. For Au 
we use the value Wco = 5.15 eV recommended by 
Michaelson [2^. The Au^~ clusters appear to be stable 
when TV > 20. However, with another value Wco = 4.3 
eV, proposed by Fomenko ^] , Eq. ifS^ gives Aujf and 
AuJj^Q critical clusters. Specific features in the energy 
properties of gold clusters were noted by Garron [5]| . 

Finally, we apply our computation procedure to the 
case of positively charged clusters NaJ^"*". Naher et al. 
|24| determined experimentally critical numbers N — 64, 
123, and 208 for clusters with n = 3, 4, and 5, respec- 
tively. In our model, we have N = Ni + AN* , and 
n = AA^*. The critical sizes are calculated using Eqs. 
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Figure 5: Size dependences of the ionization potential 15211 
Au^~ and Ag^~ clusters (solid line). Dashed lines represent 
the quasiclassical dependences IP{N) 1361 1. Arrows indicate 
the experimental critical numbers A*"* . 
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Figure 6: Critical sizes of positively charged clusters NaJ^^. 
Solid line - quasiclassical dependence 1)41^ . dashed line - the 
Rayleigh dependence, (■) - experiments, and (□) are the 
results of quantization. 



(BT|l - in which we replace Wco —E^'~' and elim- 
inate Hei- We use £coh = 1-13 eV and IP{1) = 5.14 eV 
in the calculations. 

Our results for the NaJ^+ clusters are presented in Fig. 
El For the small sized clusters, the situation can be 
described by the Rayleigh's formula according to which 
IQrI oc 7V^/^. The quasiclassical instability leads to the 
relationship eZAN^ cx N^^^. Calculated critical sizes of 
the clusters occur to be overrated as compared to exper- 



imental values. From the data presented in Fig. 5, one 
can assume that cubiform clusters transforms predomi- 
nantly into cuboids. The change in energy of the highest 
occupied state E^'^ is significantly less than that associ- 
ated with the charging due to an increase in capacitance. 
Elongation of clusters is accompanied also by the change 
in the cohesion energy (see Eqs. (|41|) and (|43|) '). This is 
confirmed by experimental data on deformation of point 
contacts. A decrease in the dimensiona lity gives a consid- 
erable increase in strength of a contact [IJ • These factors 
can be responsible for the difference between dependences 
AN*{N) calculated and experimentally found. 



V. SUMMARY 

The problem of energetics of finite metallic systems has 
been considered. We have developed an analytical the- 
ory of size dependent energy and force characteristics for 
metal slabs by using rectangular-box model, free-electron 
approximation and the model of a potential well of fi- 
nite depth. This approach makes it possible to calculate 
the oscillatory size dependence of the work function. It 
is concluded that the work function of low-dimensional 
metal structure is always smaller that that of the semi- 
infinite metal. The theory has been applied to calculate 
the elastic effects in thin metallic slab. 

In the framework of the two-component model, the size 
dependent Coulomb instability of charged metallic clus- 
ters has been described. A mechanism of the Coulomb 
instability in charged metallic clusters, different from 
Rayleigh's one, has been discussed. We showed that the 
two-component model of a metallic cluster in quasiclassi- 
cal approximation leads to different critical sizes depend- 
ing on a kind of charging particles (electrons or ions). For 
the cubiform clusters, the electronic spectrum quantiza- 
tion has been taken into account. The model enables us 
to discover the physical origin of the instability and to ex- 
plain the critical sizes of charged Au^~ , AgJ^~ and NaJ^"*" 
clusters of different shape. Results of this investigation 
might have an application in diagnostics of ultradispersed 
media, single-electronics, and in nanotechnology. 
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